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QUANTIFYING INSTABILITY SOURCES IN LIQUID ROCKET ENGINES 
Introduction 

An investigation to develop computational methodology to predict the effects of 
combusting flows on acoustic pressure oscillations in liquid rocket engines (LREs) was 
conducted. The methodology is intended to explain the causal physics of combustion driven 
acoustic resonances in LREs. A rigorous solution of the turbulent, multi-phase, multi- 
dimensional, transient mass, momentum, and energy conservation equations with combustion 
provides all of the information required to describe acoustic resonances in rocket engines. 
Ideally, a computational fluid dynamics (CFD) code would provide such a solution. However, 
real rocket engines are so geometrically complex and contain so many internal flow control 
devices that such solutions have not been produced. With current computer power, such 
solutions are becoming possible. This work is an initial attempt to predict the acoustic field in 
a simple, but typical, liquid rocket engine (LRE). Classically, acoustic phenomena are simulated 
by solving a Poisson equation for pressure by over simplifying the description of all other 
transport phenomena until such a solution is possible. 

The approach used in this work was to perform a CFD analysis of the combustion 
process with sufficient temporal resolution to identify the pressure/combustion interactions. A 
test case was chosen, which contained the physical elements known to be important in 
influencing stability, for developing the required boundary conditions and physics for a CFD 
solution to be useful. The commonly used impinging type injector elements are a major source 
of instabilities and were considered to be essential in the chosen test case. Since liquid 
propellants cannot practically be pre-mixed, mixing in the combustion chamber must be 
described. The nozzle throat serves as a well defined downstream boundary condition. 
However, the upstream boundary condition must be specified at a point which precludes 
upstream travel of pressure waves. It is not feasible to choose a generic boundary condition 
because there are so many variations of the engine feed and injector systems. Rather a specific 
configuration was investigated for which test data are available to verify the solution. 
Unfortunately, a typical engine experiment which was geometrically simple enough to be 
attractive for CFD analyses was not found. The best compromise between accepting some 
geometric complexity while retaining impinging injectors and observing instability phenomena 
was determined to be one of the subscale, Fastrac gas generator combustors investigated by 
MSFC. This combustor and operating conditions were used as the CFD test case. 

Analysis of a Fastrac Gas Generator Type Combustor 

Data from a gas generator combustor tested in the Fastrac development which manifested 
longitudinal instabilities was provided by Rocker [1]. The exit flow was choked; the injector 
elements were triplets. This case (P3 15-97-04, 6/10/97) was chosen for our investigation. A 
further attractive feature of this combustor (See Figure 1) was that the injector elements were 
located with fairly good circumferential symmetry. Two rows of F-O-F triplet elements were 
used, and the two row pattern was repeated circumferentially around the injector faceplate. A 
small nitrogen flow was injected along the centerline. The feed streams were controlled by 
cavitating venturi meters about 14 feet upstream of the faceplate. Hence, the flow was 
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controlled by the feed system cavitating venturis and injector element orifices and by the 
combustor nozzle throat. To expedite this brief Phase I investigation, an analysis of the feed 
system was not performed. Also, the combustor flowfield was assumed to be axisymmetric. 
The propellants were LOX and RP-1 . The operating conditions were such that the mixture ratio 
used was low enough to form only a small amount of soot. The approach used for the 
investigation was to simulate the unsteady case with low temporal resolution to obtain an 
essentially steady-state solution, then to increase the temporal resolution to reveal instability 
phenomena. The feed system flowfields were not sonic*, therefore, the possibility of injector 
flow chamber flow coupling must be considered. However, for these Phase I analyses, only 
chamber flow was simulated. 

The CFD code used in this study is the FDNS-RFV code. This code uses real fluid 
properties to represent the propellants, be they gaseous or liquid. Thus, the propellants enter 
the combustor as liquids, as they mix react and heat-up they become gases. The flowfield so 
predicted corresponds to a homogeneous spray combustion. The gases and liquids are modeled 
as being locally at the same temperature and pressure. The atomization and vaporization are 
thereby treated as being very fast processes. The mixing between the propellants is modeled 
with the turbulent viscosity submodel, which is the k-e model. The vaporization could be treated 
as a finite-rate process and the turbulence mixing submodel could be modified, if test data to 
quantify these processes were available. The combustion submodel pyrolyzes the RP-1 to a C 2 
intermediate which in turn combusts with oxygen to CO and H 2 . The combustion is completed 
with the wet CO mechanism. Soot is also allowed to form and to bum. Details of this 
combustion model are given in Table 1 and Ref. [2]. The rate-equations in the combustion 
model were tuned to give good predictions of temperature, molecular weight, and soot 
concentrations which have been observed in gas generator experiments [3, 4]. Matching these 
parameters should give accurate predictions of sound speed in typical gas generator type 
environments. 

To simulate the combustor flow as an axisymmetric flow the propellants were modeled 
as coming from narrow slits at the location of the injector faceplate orifices. The steady state 
CFD solution is obtained by specifying the inflow, estimating the remainder of the flowfield, 
and iterating until the exit flow equals the inflow. Test data for the chamber pressure and 
temperature were used to verify the solution. Since the specified initial flowfield is significantly 
different from the solution, large pressure waves are initially produced. These waves decrease 
in amplitude until the steady conditions are reached. The final phase of the solution produces 
a slow, monotonic approach to the final answer. The final solution is shown in Figure 2. The 
predicted chamber pressure was 465 psia, compared to the 470 psia measured value. The 
measured flowrate was 7.05 lb/s; that from the simulation was 6.87 lb/s. The predicted 
temperatures are very close to the measured values of 1545-1573°R. The predicted mean flow 
properties near the end of the chamber are shown in Table 2. For the predicted sound speed in 
the chamber and the chamber geometry, the first longitudinal mode would correspond to 798 Hz. 
Equilibrium chemistry, 1-dimensional solutions were also run for this test case. Two chemistry 
models were used, one with soot formation and one without. The case with soot formation 
indicated that 44 mole% soot was formed, the chamber sound speed was 2358 fps, and the 
corresponding 1L frequency was 1278 Hz, molecular weight was 26, the temperature was 2023 
R, and the flowrate was 4.5 lb/s. The case without soot gave the chamber sound speed of 1571 
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Table la. Combustion Kinetic Models for RP-l/0 2 



Rate = AT B exp{-E/RT} [C 2 H 4 ] 0 5 [0 2 ] 10 

A = 1.29xl0 15 , B = 1, E/R = 2.516X10 4 


48 C 2 H 4 -> Soot + 84 H 2 

Rate = AT B exp{-E/RT} [C 2 H 4 ] 2 0 [0 2 + e] 0 5 
A = 5.1308xl0 12 , B = -2, 

E/R = 1.61 lxlO 4 , e = 0.001p/Mo 2 g-mole/cm 3 

Soot is defined as this approximation makes incipient soot particles the right size and 

maintains the right carbon-to-hydrogen ratio for soot. 


Table lb. Combustion Kinetic Model for Wet CO 
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Table 2. Flow Properties at the End of the Combustor from the Two-Phase 
Simulation Which Were Used to Define the Ideal Gas 


Temperature 

(°R) 

Molecular weight 
(lb m /lb mole ) 

7 

Density 

(lb m /ft 3 ) 

Speed of Sound 
(ft/s) 

1680 

41.97 

1.09 

1.106 

1472.7 

| Species Concentrations (mass fractions) 


h 2 o 

CO 

C0 2 

h 2 

c 2 h 4 

RP-1 

Soot 

n 2 

0.063 

0.2087 

0.0419 

0.0121 

0.0705 

0.5869 

0.0162 

0.0007 


fps, the corresponding 1L frequency was 852 Hz, molecular weight was 21, the temperature was 
1134 R, and the flowrate was 6.9 lb/s. The observed IL frequencies varied from 789 to 821 Hz. 
The assumption of equilibrium soot production is obviously poor. The CEC solution without 
soot production provides a reasonable estimate of the mass flowrate, but the predicted 
temperature and molecular weight differences between the CEC equilibrium and CFD solutions 
are significant. The equilibrium analysis pyrolyzes all of the hydrocarbons to small species. 
This is an endothermic process which cools the flow. The CFD analysis leaves much of the 
unreacted hydrocarbons present as vaporized RP-1. The presence of the RP-1 is necessary to 
match the mean molecular weights reported for gas generator tests [3]. Therefore, the finite-rate 
combustion model used in the CFD code is believed to offer the best estimate of flow properties 
in the Fastrack gas generator. Furthermore, the properties shown in Table 2 provides a good 
estimate of an ideal gas surrogate to represent the combusting propellants. 

Tr ansi ent Engine Simulations - Real Fluids 

The plan was to continue the steady-state analysis with the code tuned to represent time 
accurate solutions using the real fluids option of the CFD code. First, the unsteady flowfield 
solution behavior was investigated by reducing the time-step size and using the implicit time 
integration option and continuing the calculation from the steady-state solution. Figure 3 shows 
the results of this simulation. As expected when the smaller time step and the implicit 
integration scheme are used, a perturbation resulted in the flow. This perturbation is seen to 
damp out and produce another steady-state answer which differs slightly from the solution shown 
in Figure 2. Such behavior suggests that the observed instability was due to coupling with the 
feed system and/or injector element orifice flows. Since the real fluids simulation is 
computationally intensive, it was deemed impractical to extend the computational domain 
upstream to the cavitating venturi meters in this short Phase I investigation. Rather, the inlet 
flow rate was perturbed in a sinusoidal manner to determine how the pressure field would 
develop. This forced perturbation crudely represents what would happen due to feed 
system/injector created flow disturbances. Figure 4 shows that the pressure (at the impingement 
point, #1) initially responds out of phase with the flowrate oscillation, and eventually becomes 
in phase with flowrate after several cycles of oscillation. Pressure monitoring station #3 is at 
the pressure measurement station which is about 3/4 of the way down the barrel section of the 
combustor, while station #2 is located between the turbulence ring and station #3. The pressure 
oscillations at the 2nd and 3rd monitoring stations become out of phase with those at the inlet 
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after the oscillations become regular. All of the pressure oscillations tend to become small in 
amplitude, i.e. about 5 psia. To determine if the pressure oscillations would continue after the 
flowrate/pressure at the inlet boundary was set to an empirical correlation, the simulation was 
continued with this boundary condition. This correlation will be described in the next section. 
Figure 5 shows that the simulation tends to diverge, but it took an excessive amount of clock 
time to make this determination. In fact so much time, that the correlation could not be tuned 
to obtain reasonable parameters in the correlation. Two other factors could also be causing 
pressure oscillation damping: (1) the simulation could be over damped by the turbulence model 
and/or the artificial (numerical) damping in the solver, and (2) the boundary conditions used to 
in the chamber. Such factors could only be evaluated parametrically. However, again due to 
the short time frame of the Phase I investigation, the real fluids CFD code required too much 
time to provide parametric studies of these effects. Therefore, these effects were studied with 
an ideal gas model to represent the chamber flow. 

Transient Engine Simulations - Ideal Gases 

In order to determine the cause of the pressure oscillation damping in the two-phase, real 
fluid combustion flow simulation, a series of transient, ideal-gas simulations was conducted. 
The propellants entering the combustion chamber were assumed to be premixed and pre- 
combusted and have the same flow properties as the mean flow at the end of the combustion 
chamber in the two-phase flow calculation. These flow properties are listed in Table 2. As 
previously noted, the temperature, mean molecular weight, and sound speed of the pre- 
combusted gas has a nature frequency which is close to the 1-L mode of the pressure oscillation 
in the Fastrac Gas Generator. Observing from Figure 2 that the hot-gas exit flow is 
representative of a large part of the combustor flow, the ideal-gas simulations should have 
approximately the same residence time in the combustor as the real gases. 

The purpose of using the ideal-gas flow simulation to study the instability is to 
substantially reduce the computation time, such that various boundary conditions and numerical 
schemes can be tested in a timely fashion. In addition, the premixed, hot-gas flow simulation 
excludes the effects of combustion instability and spray/atomization/mixing instability, and thus 
the chamber pressure oscillations are due to acoustic waves only. 

A steady-state calculation with a uniform profile of flow properties at the chamber inlet 
was performed first to establish a converged flow solution for the chamber and nozzle. The time 
history of the chamber pressure is plotted as shown in Figure 6. In the next step, the inlet hot- 
gas mass flow rate was perturbed according to the following function: 

m = rho [ 1 + <p m sin ( 2nc^,t ) ] 


where <f> m = 0.3 and o> m = 1/800 sec' 1 were tuned to establish a 35 psi (Ap/p c * 0.074) pressure 
oscillation at the inlet with a frequency of 800 Hz. Numerical calculations with perturbed inlet 
mass flow rates were conducted for more than 10 cycles. The numerical results indicated that 
the magnitude of pressure oscillations near the chamber exit were smaller than those at the inlet. 
To insure that the attenuation of the pressure oscillation was not caused by numerical damping 
and the over-prediction of the eddy viscosity, this ideal-gas simulation was repeated with: (1) 
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lower inlet turbulence intensity, (2) the temperature-corrected turbulence model, and (3) the 
exclusion of artificial damping terms in the FDNS-RFV code. The numerical results still 
indicated that both the mean pressure and the magnitude of pressure fluctuations vary along the 
streamwise locations, as shown in Figure 7. The location of Pressure #\ corresponds to the 
impingement point in the two-phase, spray calculation. The location of the pressure 
measurement near the chamber exit is identified as Pressure #3; while Pressure #2 is located 
between the turbulence ring and the location of Pressure #3. The predicted mean pressure is 
consistent with flow physics: highest pressure and lowest mean velocity at the chamber inlet; 
lowest pressure and highest mean velocity at Pressure #2 where recirculation zone behind the 
turbulence ring reduces the mean flow area; lower pressure and higher mean velocity at Pressure 
#3 which is near the end of the flow recirculation zone. These velocity and pressure fields are 
shown in Figure 8. It is to be noted that the phase lag between different pressure points remains 
constant, which indicates the implicit time integration scheme employed by the flow solver is 
accurate. Furthermore, the pressure at the impingement point, which is very close to the inlet, 
oscillates almost in phase with the inlet mass flow rate. Since the predicted eddy viscosity is 
believed to be small enough (the reattachment length behind the turbulence ring is longer than 
8 step heights of the turbulence ring), the loss due to the viscous force should be small. Hence, 
the reduction of the chamber pressure and its fluctuations should be real phenomena. The most 
telling evidence that the variation in the chamber pressure amplitudes with location is real is that 
the smallest amplitude is at station #2. If the amplitude variation were due to dissipation only, 
the amplitudes could not decease to a minimum and then increase along the combustor. These 
observations suggest that relatively minor geometric complexities can have a profound influence 
on pressure waves in liquid rocket engines. 

Several transient calculations were performed by starting with the predicted perturbed 
flowfield and by not using the inlet mass flowrate perturbation. Since the numerical calculations 
started at the injector face, the effects of upstream pressure/mass flux perturbations from either 
the injector orifices or the feed lines were not simulated. Numerical experiments were 
conducted to find a set of inlet boundary conditions, with which the chamber pressure 
fluctuations dissipated the least and the mean chamber pressure at specific locations remained 
constant during the transient numerical calculations. The inlet boundary condition for the static 
pressure is either zero gradient or constant gradient. There are four types of inlet boundary 
conditions for the velocity: constant mass flow rate, constant total pressure, constant velocity, 
and a correlation of the inlet mass flowrate fluctuation to the inlet pressure fluctuation which can 
be expressed as, 

Am Ap 

m p P 0 

(f> is a constant tuned to meet the goal mentioned above. It was found that the constant pressure 
gradient boundary condition combined with any of the other velocity boundary conditions tends 
to overshoot or undershoot the mean chamber pressure. For example, see Figure 9. In 
addition, the pressure oscillations at various locations were damped out with time. The 
combination of zero pressure gradient with the correlation of mass flux fluctuation to the 
pressure fluctuation showed a larger variation in the numerical results for different values of <f> p . 
For negative <f> p (-5 and -1) where the oscillations of inlet mass flow rate are opposite to the 
oscillations of the inlet pressure, not only the mean pressure is lower but also the pressure 


6 



SECA-FR-00-04 


fluctuations decay faster compared to the results of positive <f> p , as shown in Figures 10 and 11, 
respectively. Moreover, the pressure fluctuations predicted with 4> p =-5 decay faster than those 
predicted with <f> p —-\. For positive 4> p (5 and 1) where the oscillations of inlet mass flow rate 
are in phase with the inlet pressure oscillations, the pressure fluctuations also decay but with a 
slower rate, see Figures 12 and 13. However, the mean pressures at various locations diverge 
if 4> p =5; but the mean values remain constant if <£ p = 1. The combination of zero pressure 
gradient with constant inlet total pressure produced a high frequency mass flux fluctuation, still 
the pressure oscillation decayed very fast, as shown in Figure 14. The zero pressure gradient 
boundary condition in combination with either constant inlet velocity or constant inlet mass flux 
displayed the same characteristics as can be seen in Figures 15 and 16, respectively. Mean 
pressure at a given location remains constant with time, and the pressure fluctuations decreases 
with time. It should be noted that the decay of pressure fluctuations predicted by the 
combination of zero pressure gradient and constant inlet mass flux is the slowest compared to 
the other sets of inlet boundary conditions. 

Based on these simulations, it is concluded that the combustion chamber and the nozzle 
will not be able to sustain instabilities. Since the pressure perturbations decrease along the 
streamwise direction due to the viscous and turbulence effects and the complex flowfield which 
exists in the chamber, smaller pressure waves return to the inlet boundary from the nozzle, and 
subsequently these pressure waves damp out. By examining the pressure test data, the pressure 
waves seem to exhibit a behavior such that the pressure oscillations decay for 3 cycles and are 
followed by a pressure spike. This characteristic is repeatable throughout the test data. It is 
suspected that there is a long-wave perturbation from upstream (feed line and/or injector orifices) 
which maintains the instability. It is further concluded that the observed instabilities in this 
Fastrack gas generator configuration are due to coupling between the feed line and/or injector 
element orifice flow coupling and the combustor flow proper. Further computational analyses 
should include a domain which extends through the injector orifice and to the cavitating venturis. 
The combustion physics included in the two-phase, real-fluid combustion models in the FDNS- 
RFV code should be used for these simulations even though they are time consuming to perform. 

Conclusions and Recommendations 

1. Real fluid CFD simulations of combustion stability phenomena in LREs are very 
computationally intensive, and should be addressed with the fastest processors available. 

2. The rapid flowfield damping of real and ideal gas simulations of Fastrack type gas 
generators indicates that there is strong coupling between the feed line/injector element 
flowfields and those flowfields in the combustion chamber. 

3. Unnecessary artificial damping in the FDNS-RFV code has been eliminated to 
expedite combustion stability simulations. 

4. Numerous methods of applying chamber boundary conditions which were investigated 
did not significantly affect the dissipation of the perturbed solution to a steady or no-flow 
condition. Some of the boundary condition specifications considered caused an immediate 
computational instability to occur. 
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5. Oscillatory flows driven by chamber boundary condition specifications were shown 
to be modulated and for pressure waves to form because of local flowfield variations. Such 
variations are due to configuration variations as well as flow dissipation. 

6. Future CFD stability analyses of Fastrack type gas generator combustors should be 
conducted for a computation region which includes flow control devices on both the in-flow and 
out-flow boundaries of the region. In-flow controls may be of the cavitating venturi type and 
out-flow, sonic throats of the nozzle. 

7. Parallel processing on multiple PC processors should be used to provide the computer 
power needed to more completely address combustion stability issues of liquid rocket engines. 
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Flow Properties of the Fastrac Gas Generator 




Pressure Oscillations of the Fastrac Gas Generator 

(Real Fluid, Non-Perturbed) 
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Figure 14 

Pressure and Mass Flow Rate Oscillations of the Fastrac Gas Generator 
(Idea Gas, Perturbation Released, Constant Pressure Gradient, Constant Total Pressure) 
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Figure 16 

Pressure and Mass Flow Rate Oscillations of the Fastrac Gas Generator 

(Idea Gas, Perturbation Released, Constant Pressure Gradient, Constant Mass Flux) 
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